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I. INTRODUCTION 



Transport properties of deterministic chaotic Hamiltonian dynamics can exhibit a rich 
variety of phenomena. In many situations where chaos is strong, transport through phase 
space tends to be uniform and can imitate random processes of the Brownian type. The 
motion on average is a normal diffusion. However, in more generic situations where regular 
as well as chaotic behaviour coexist the phase space has a fractal structure, and this structure 
can induce important variations of the transport process. 

In this paper our purpose is to investigate transport properties of a deterministic two- 
dimensional area preserving map whose dynamics can be either mixed or strongly chaotic 
as a parameter is varied. It is defined by the time- dependent Hamiltonian 

H(q,p,t) = - V 2 cos(27rp) - Vi cos(27rg) K{t) 
K(t) = TY^-ooKt-nr) 
called the kicked Harper model. It arises in several branches of physics, as for example in the 
description of charged particles interacting with an electrostatic wave packet in a transverse 
uniform magnetic field fl|-|3|], in the study of the two-dimensional motion of electrons in a 
periodic potential plus uniform perpendicular magnetic field [§J and in connection to the 
steady-state dynamics of two and three-dimensional ideal incompressible fluids A 
characteristic feature of Harper-like Hamiltonians is the particular functional form of the 
kinetic energy. The equations of motion induced by Eq.([|) can be integrated over one period 
r of the kicking potential, leading to the map 

p n+1 = p n - 71 sin(27rg n ) 

(2) 

q n +\ = q n + l2 sin(27rp n+ i) 
where 7$ = l-nViT. We can consider the motion either in the entire phase-space plane or in 
the fundamental domain (q,p) G [0, 1[ 2 . The existence of a fundamental domain comes from 
the fact that the map (|2]) is invariant under translations (q,p) —>■ (q + l,p + m) - where (I, m) 
are two arbitrary positive or negative integers - since such a shift of the initial conditions 
produces a global shift of the trajectory. 
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In the limit r — > the model reduces to the continuous time evolution of the Harper 
Hamiltonian H (q,p) = —V% cos(27rp) — V\ cos(27rg) which, from the dynamical point of view, 
is an integrable system. For simplicity, we restrict the following discussion to the special case 
7i = 72 = 7 (V\ = V2). Then the phase space plots of H are characterized by the existence 
of a separatrix defined by H (q,p) = connecting both sides of the fundamental domain. 
In the full phase plane, the separatrix forms a regular lattice having a square symmetry (cf 
Fig. 1(a)), while the contour lines Ho(q,p) = const 7^ are closed curves. For 7 > the 
separatrix is destroyed and becomes a mesh of finite thickness inside which the dynamics 
is chaotic. This chaotic mesh, which has been called the stochastic web QTJ] , coexists with 
regions of phase space where regular motion takes place. As 7 increases the stochastic web 
grows, and for sufficiently large values of 7 the web disappears and chaotic dynamics covers 
most of phase space. This transition is shown in Fig.l for the motion restricted to the 
fundamental domain. In order to better visualize the geometry of the stochastic web, in 
Fig. 2 we have plotted again Fig. 1(b) but now in the full phase plane. 

The motion of the particle is not chaotic for any initial condition inside the stochastic 
web, since there are regular islands on it (clearly visible in Fig. 2). However, if an initial 
condition corresponding to a chaotic motion on the web is chosen then, because the web 
extends over the whole phase-space plane, the particle can "diffuse" arbitrarily far from its 
initial position. In the p direction, the diffusion will be characterized by the coefficient 

D = lim ^-P,?) R (3) 

n— >oo fi 

where the average is taken over a set R of phase-space initial conditions. This coefficient can 
in general be expressed as an infinite series involving the force-force correlation functions (see 
and Eq.([8]) below). The existence of D (i.e., the convergence of the series) and therefore 
the existence of a true deterministic diffusive process is in general difficult to prove rigorously. 
This, however, has been done for a different map, the sawtooth map J7|. Moreover, several 
aspects of the diffusion process in deterministic maps have been studied in the past by 
truncation of the series and/or adding some noise. These studies were mainly concentrated 



on the standard map [pHIOfl or on different types of stochastic webs JT|,[TT|. In particular, 



computations of D for the kicked Harper map can be found in fi"2| , |T3|| . Our aim here is to 
compute the diffusion coefficient as a function of 7 up to second order in the correlation 
functions and to analyze its dependence on 7 in terms of the underlying classical dynamics 
(section III). We obtain, for all 7, a finite value of D which is in very good agreement with 
the numerical results except near some parameter ranges where accelerator modes appear. 
These modes are stable periodic orbits which coherently propagate across the web. The 
shortest periodic orbits corresponding to accelerator modes as well as their stability are 
computed in Section II. The existence of these orbits produces a divergence of the diffusion 
coefficient if the average in Eq. @ is taken over all possible initial conditions. This is because 
in the presence of accelerator modes the average can be splitted into two parts: the first 
is over the stable phase-space regions <S where accelerator modes exist, while the second is 
over the chaotic component C 

((Pn - Po) 2 ) Rf = ((Pn- Po) 2 )s + ((Pn ~ Po) 2 )c • (4) 

where the index Rf means that the average has been taken over the whole fundamental 
domain. Because accelerator modes propagate in a ballistic manner, ((p n —po) 2 )s ~ n 2 and 
the diffusion coefficient associated to this region diverges like n. 

The presence of accelerator modes can have more unexpected consequences on the trans- 
port properties of the system than the one just mentioned. Their influence extends into the 
chaotic sea and may produce non-local effects. In fact, the partial barriers (cantori) induced 
by their presence in the chaotic regions of the web may modify the nature of the transport 
in the stochastic web as well. For most of the time the particle diffuses, but sometimes it is 
trapped in the neighbourhood of an accelerator modes and drag along in a ballistic manner. 
These rare events can dominate the dynamics and produce anomalous diffusion of the Levy 



type P 



((Pn-Po) 2 )c~^, l</3<2. (5) 
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As a consequence the diffusion coefficient associated to the chaotic component also has a 
power-law divergence with exponent (3—1. 

Anomalous diffusion was already observed numerically in different Hamiltonian systems, 
like the standard map p5|-p^], in stochastic webs [0,23] and in two-dimensional smooth 
potentials [j21{| . There are however few systematic quantitative results. At present there is 
no theory making the connexion between the existence of phase space barriers in the chaotic 
layers on the one hand and anomalous diffusion of the Levy type on the other hand. Only 



in one case it was possible to establish [p2| a direct connexion between some parameter of a 
map and the temporal exponent (3 determining the anomalous diffusion. But unfortunately 
this was done in the simple case of a one- dimensional iterative map. In Sect. Ill we compute, 
as a function of 7, the exponent (3. Our numerical results, presented in Fig.5 and Fig.6, 
demonstrate that the kicked Harper map provides an additional example of anomalous 
diffusion of the Levy type in Hamiltonian systems. Different aspects of this exponent in 
connexion to the existence of normal modes are discussed. 



II. ACCELERATOR MODES 

The period-1 fixed points of the map @j are determined by the condition p n+ \ = p n , 
Qn+i = Qn- We are interested however in the fixed points of the map on the fundamental 
domain, and thus we identify phase-space points which differ from an integer value: p„+i = 
p n + w p , q n+ i = q n + w q leading to 

sin(27rp ) = -w q /ji 

(6) 

sin(27rg ) = -w p /j 2 

where (w q , w p ) 6 Z 2 and we have used the notation q n = q , p n = p . These two integers are 
the winding numbers of the period-1 fixed points since they define the number of turns the 
orbit does around the fundamental domain in one iteration of the map. In the phase plane, 
periodic orbits with winding numbers (w q , w p ) propagate ballistically, with (p n — p ) = w p n 
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and (q n — go) = w q n. But not all of them correspond to accelerator modes either because 
they are not stable or because they do not propagate (w g — w p — 0). 

For simplicity, we again restrict the discussion to 71 = 72 = 7. For a fixed value of 7, the 
solutions of Eqs.(||) are characterized by all integers such that \wi\ < [7], where the square 
brackets mean integer part For each possible set (w q , w p ) of winding numbers there are 
four different solutions in [0,1[ 2 : (<7o,Po), (go,l/2-po), (V 2_ 9o,Po) and (1/2- q , l/2-p ), 
where go = arcsin(— w p /^)/2ix and po = arcsin(— w q /j)/2n. For < 7 < 1 only (w q ,w p ) = 
(0, 0) is a solution and as a consequence there are no period-1 accelerator modes in this 
parameter range since these orbits do no propagate. As 7 increases, at each integer value 
there is a first-order bifurcation (fold catastrophe) where several new orbits appear. Taking 
into account all possible new couples of winding numbers and their signs (where at least one 
is equal to [7] while the other one can be less than [7]), there are at each integer value of 
7 a total number of 32 [7] new orbits appearing. This makes, for a given value of 7, a total 
number of (4 [7] + 2) 2 period-1 fixed points. 

These orbits are stable if |TVM| < 2 and unstable otherwise, where M is the monodromy 
matrix. For a given couple (w q ,w p ) the trace is given by 

TrM = 2 - 4vr% pV /( 7 2 - % 2 )(7 2 -^) (7) 

with e qp = sgn(cos(27rg ) cos(27rp ))- At any integer value of 7, TrM = 2 for all the periodic 
orbits emerging at that value of the parameter. e qp is positive for half of them, and negative 
for the other half. For the latter half, it follows from fl7|) that TrM > 2 as 7 increases from 
the bifurcation point. On the contrary, those having e qp = 1 are stable as 7 increases from 
[7] but as soon as TrM < —2 they become unstable too. This range of stability defines 
the domain of existence in parameter space of period-1 accelerator modes. Its width is in 
fact relatively small (compare to one), and its size decreases with [7]. For example in the 
particular case w q = w p = w = [7], TrM = -2 at 7 C = Q7] 2 + 1/tt 2 ) 1/2 w [7] + l/(2vr 2 [7]). 

At the value 7 C where TrM = —2 there is generically a period doubling bifurcation, 
where a period-2 stable fixed point of the map appears. If 7 further increases, this or- 
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bit will soon become unstable through another bifurcation giving rise to the usual cas- 
cade of bifurcations. We are not going to establish here a complete classification of the 
period-2 fixed points, but briefly discuss a special family who plays an important role 
in the variations of the diffusion coefficient. These are the fixed points that emerge at 
half integer values, 7 = m/2, m = 1,3, Their point of emergence is independent 
of 7 and given by (?o,Po) — (1/4, 0), (3/4, 0), (0, 1/4) and (0,3/4) with winding numbers 
(w q , Wp) = (0, —m), (0, m), (m, 0) and (— m, 0), respectively. As for the period-1 orbits, their 
position and stability varies with 7 and there is a small range in parameter space close to 
7 = m/2 in which some of them are stable. Accelerator modes of period 2 therefore exist 
around all half integer values of 7. The analysis of the accelerator modes can be pursued 
for arbitrary high periods, with increasing complexity. 



III. DETERMINISTIC DIFFUSION 



From Eqs. (0) it follows that 



n-l 



Pn = Po - 7i sin(27r<?i 



i=0 



and then 



(Ap n ) 2 = (p n -p ) 2 = 7i 2 



n— 1 n—2 n—l 

^2, sin 2 (27rgj) + 2 ^ ^ sin(27rgj) sin(27rgj 

i=0 j=i+l 



i=0 



Our purpose is to compute the diffusion coefficient (|) averaged over all possible initial 
conditions (qo,Po) in the fundamental domain. By definition 



</>*, 



1 rl 



JO 



dq dp f{qo,Po) ■ 



It is therefore implicit when computing any average value that the coordinates and momenta 
(qi,Pi) are expressed, using the equations of motion (g), as a function of (qo,Po)- Because 
(sin 2 (27rgj)) = 1/2 and introducing the force-force correlation function 



C(t) = (sm(27iqi) sm(27iq i+T )) Rf 



it follows that 



n-1 



A 



n 



T = l T = l 

where we have used the invariance of C(r) under time translations. If the C(r) fall off 



rapidly enough with r, then |7[] 



g=lim «W>i,, = 2? 

2 



l+4£C(r) 



T = l 



(8) 



If correlations are neglected, i.e. the are considered as random uncorrelated variables, 
then it follows from Eq.(Q) that the motion in the p-direction reduces to a random walk 
with length-step "f/y/2. Accordingly, we obtain from Eq.(H) D = "fx/2. This is expected 
to be a good approximation for large values of the 7's where the "randomization" is more 
effective. From this point of view, the corrections for C(r), r = 1,2, ■ • ■ in Eq.(§) can be 
considered large-7 expansion. 

The following relations are needed in order to compute the correlation functions 



(9) 



cos(zsm9) = Jp(z) + 2 JikiA cos(2k9) 

k=l 

00 

sin(zsin^) = 2 ^ J 2k +i{z) sin((2A; + 1)0) . 

k=0 

The one-step correlator is 

C(l) = (sin(2vrg ) sin(2vrgi)) R/ 

= (sin(27rg ) sin [2iTq + 2^72 sin(27rp - 27r7i sm(2nq ))}) R f 

From Eqs.(^) and simple trigonometric relations the averages can be computed explicitly. 
One finds 



C(l) = ij (27r 72 ) 



(10) 



The computation of the correlator for points differing by two time steps 

C{2) = (sin(27T(Z_i) sin(27r 9l )) % 

= (sin [27rg — 27T7 2 sin(27rp )] sin [27rg + 27T7 2 sin(27rp — 2^7! sin(27rg ))])_R / 
follows the same lines, with increasing algebraic complexity. We get 
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1 oo 

C(2) = - J 2 (27r 72 ) - £ [(-I) n+ V (2n7r 7 i) + J 2 (2mr 7 i)] J 2 (2tt 72 ) . (11) 

Z n=l 

From the computation of C(3) we only keep the term Jq (27r 7 2)/2, the other being of higher 
order in powers of Bessel functions. Then, from Eqs. ( fTO] ) and ( |TTD we finally have 

D/^l/2) « 1 + 2[J (27T 72 ) + J 2 (2tt 72 ) + J 3 (2tt 72 )] 

( 12 ) 

-4£ [(-l)" +1 Jo(2n7r 7l ) + J 2 (2vm 7l )] J 2 (2vr 72 ) . 

n=l 

In order to compare this result with the actual properties of the transport, we have numeri- 
cally computed the quantity _D/( 7 2 /2) as a function of 7 = 7 i = 72 . For each value of 7 , the 
diffusion coefficient was calculated from the slope of ^^ 2 ) f° r n = 0> * ' * > 3000 averaged 
over 9000 different initial conditions taken at random on the elementary cell. The values 
obtained are plotted in Fig.3, together with the analytical result fll2"D . The approximation 
(|12|) clearly accounts for the main features of the diffusion process. Indeed, for most values 
of 7, the numerical points are in very good agreement with the finite value obtained in 



(|12|). The oscillating behaviour of period A7 = 1 can be physically related to the strong 
modifications induced in the transport properties by the appearance of period- 1 accelerator 
modes at each integer value of 7. For large values of 7 the oscillations are more regular with 
decreasing amplitude oc A/tt^/^. For small values of 7 the phase space structure is more 
rich and accordingly D is less regular. In this regime the agreement with Eq. (|T2|) is less 
good, as shown in the inset of Fig.3, and below 7 rs 0.35 the theoretical value is definitely 
wrong (indeed, from Eq.([l2|) we get D / /2)\ 1= q = 7, while the numerical values tend to 



0). To compute more accurately the diffusion coefficient in this limit, we make use of the 



results obtained in []T[. In the limit 7 < 1 the diffusion coefficient of the stochastic web 
takes the form 

27 s 

Dstoch ~ tt + 27 [3-ln(c7r) + 31n(7r 7 )] ' 7 <<= 1 (13) 
where c ~ 20. Because our computations are averaged over all possible initial conditions, 
we shall multiply Eq. ( TT5| ) by the relative area A of the stochastic web jjTT] 



A « ^e~^ 7 « 1 
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with the result 



D 



4c 



(14) 



(7 2 /2) 



7r 2 7 2 {it + 27 [3 — ln(c7r) + 3 111(777)] } 



7 < 1 • 



This function is plotted at the origin in the inset of Fig.3, and is in very good agreement 
with the numerical values. 

On the other hand, the numerical results do not agree with the theory close to low integer 
and half integer values of 7, where some of the numerical points lie systematically higher. If 
the number of iterations employed in order to compute the slope is increased, these points are 
not stable but rather move up, a fact clearly pointing out towards a divergence of the diffusion 
coefficient. Analytically, the divergence should come out from further correlation terms in 
the expansion (|8]) As mentioned in the introduction, this divergent behaviour is closely 
related to the appearance of accelerator modes. Because of their stability, these ballistic 
modes occupy a (small but) finite phase-space volume. It then follows that ((p n —po) 2 )s ~ n 2 
and D — ► 00 as n when n — > 00 since the average was computed over the whole fundamental 
domain (cf Eq.(|4])). For large integer or half integer values of 7 the divergence is not visible 
numerically because the area and range of existence of the accelerator modes is very small 
(cf Section II). 

One could think that the divergence of the transport coefficient is limited to stable regions 
of phase space where accelerator modes exist. This in fact is not true. When this modes 
are present, two rather different dynamics coexist. On the one hand, there is the ballistic 
motion associated with the small stable regions. On the other hand there is the chaotic sea 
(occupying most of phase-space as soon as 7 > 0.5), where incoherent motion takes place 
and where we may expect (as in the absence of accelerator modes) a regular Brownian- 
type diffusion. But the frontier between both types of motion is not sharp. Stable islands 
are surrounded by broken tori (called cantori). These broken tori play the role of partial 
barriers for the motion in the chaotic sea. When the particle is trapped by a cantorus in 
the neighborhood of the accelerator mode, it is dragged across the phase plane ballistically. 
These rare events can however dominate the dynamics and produce an anomalous diffusion 
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23j,18fl. If the trapping time probability distribution has the asymptotic behaviour 

^(n) « l/n 4 ~P 



then we can classify the dynamics into three distinct regimes Jl8| 

n 2 2 < (3 < 3 ballistic 
((A Pn ) 2 ) ~{ n Pl<P<2 Levy (15) 
n (3 < 1 Brownian 
The presence of accelerator modes may then induce strong modifications into the transport 
properties of the chaotic component. This happens when the exponent j3 lies between 1 
and 2. The anomalous trajectories consists of a regular diffusion process punctuated by 
intermittent ballistic episodes. We find that this type of transport occurs in the kicked 
Harper map in certain parameter ranges. As an example, in Fig.4 we have plotted 5000 
iterations of a trajectory lying in the chaotic component obtained for 7 = 1.011. The 
orbit is shown in the entire phase plane as well as in the fundamental domain. The former 
picture stresses the ballistic Levy-type intermittent behaviour of the orbit, while the latter 
makes clear the trapping of the orbit around the accelerator modes. In fact, the long step 
propagating coherently over more than 150 fundamental domains is clearly visible in part 
(b) of that figure as a concentration of points around the stable region associated with the 
period-1 orbit emerging at 7 = 1 at (q,p) = (1/4,3/4). 

In order to make more quantitative these considerations, we have numerically computed 
the exponent (3 as a function of 7 for the chaotic component of phase space. This was 
done propagating for each value of 7 9000 trajectories iterated 5000 times. The average 
over the trajectories was taken in order to compute from Eq.(|5]) the exponent (3. The initial 
conditions of these trajectories were taken at random in a small region of size 0.05 around the 
phase-space point (0, 1/2). This choice ensures that the orbits lie in the chaotic component 
of phase space. The results obtained are plotted in Fig.5 for 0.2 < 7 < 2, where it is clearly 
visible that values of (3 greater than one are obtained for values of 7 slightly higher than 0.5 
and 1. This confirms the existence of anomalous diffusion. However, it is remarkable that 

11 



in spite of the fully developed fractal structure of phase space existing for < 7 < 1, we 
observe anomalous diffusion only in a very localized region around 7 ~ 0.5. Could be that 
looking with a smaller step in 7 other peaks will emerge. In Fig. 6 the two peaks of Fig.5 are 
amplified. The function ^(7) in the second peak is a strongly oscillating function, and for 
7 pa 1.05 the anomalous diffusion disappears. This value is in good agreement with the value 
7 C pa 1.051 obtained in Section II for the disappearance of the (last) period-1 accelerator 
mode emerging at 7 = 1. 

IV. CONCLUSION 

The normal diffusion process of a deterministic two-dimensional area-preserving map can 
thus be very well described in terms of a truncated expansion of the diffusion coefficient. 
The theoretical prediction works well even in parameter ranges where the fractal phase-space 
structure of the classical dynamics is strong (mixed regime). However, when stable periodic 
orbits propagating coherently through the chaotic sea are present, the diffusion coefficient 
diverges. This divergence is trivial in the sense that it is directly related to the existence of 
a finite and stable ballistic region of phase space. It is not accounted for by the truncated 
expansion but should be recovered if additional terms are added 0. 

The presence of stable coherent orbits may also modify the transport properties of the 
chaotic sea itself. We find that in certain parameter regions surrounding the emergence of 
accelerator modes the diffusion is anomalous, with 1 < (3 < 2 in Eq.fl^) (and fllBD ). The 
puzzling fact is that this is not always true. For example, we were not able to detect any 
anomalous diffusion for parameter values close to (and slightly higher than) 7 = 1.5, 2.5, • • •, 
where we know that period-2 accelerator modes emerge, but we do observed anomalous 
diffusion at 7 = 0.5, 1, 2, 3, • ■ -. This strengthens the need to develop a quantitative and 
predictive theory for the anomalous diffusion in Hamiltonian systems, which does not exist 
today. 
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FIGURES 



FIG. 1: The kicked Harper map plotted in the fundamental domain, (a) 7 = 0.00063, (b) 
7 = 0.31, (c) 7 = 0.565 and (d) 7 = 0.63. 

FIG. 2: Same as in Fig. 1(b) but plotted in the full phase plane. 

FIG. 3: Normalized diffusion coefficient as a function of 7. Analytical curve (solid line); 
numerical results (dots). Inset: amplification in the range < 7 < 2; the full curve at 
the origin represents Eq . (]T4[) . 

FIG. 4: Chaotic trajectory obtained for 7 = 1.011 where anomalous diffusion occurs, 
(a) phase plane representation; (b) same trajectory as in part (a) but plotted in the 
fundamental domain. 

FIG. 5: The index (3 of Eq.(|]) as a function of 7 (dots). The interpolating full line is to 
guide the eye. 

FIG. 6: Amplification of the two peaks of Fig.5. (a) Peak around 7 = 0.5. (b) Peak around 
7 = 1. 
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